library(dplyr)
library(lattice)
library(moments)
library(tidyr)
library(GGally)
set.seed(2)

1. Выбор данных

https://www.kaggle.com/shivam2503/diamonds

2. Считывание и просмотр

df <- read.csv("diamonds.csv", header = TRUE, as.is = FALSE)
df <- sample_n(df, 1000)
head(df)
##       X carat       cut color clarity depth table price    x    y    z
## 1 46031  0.57     Ideal     G     VS2  61.6    57  1728 5.36 5.32 3.29
## 2 11014  1.01     Ideal     D     SI2  61.3    54  4916 6.47 6.52 3.99
## 3 36044  0.45      Fair     F     VS2  67.0    56   923 4.77 4.70 3.17
## 4 15657  1.04   Premium     H    VVS2  59.1    60  6278 6.66 6.60 3.92
## 5 11851  0.90 Very Good     G    VVS2  59.8    60  5102 6.23 6.28 3.74
## 6 14800  1.20 Very Good     I     VS2  62.3    56  5955 6.76 6.81 4.23

3. Описание данных

Этот классический набор данных содержит цены и другие атрибуты почти 54 000 бриллиантов.

X— индекс

carat— вес бриллианта в каратах

cut— качество огранки

color— цвет бриллианта

clarity— включения в бриллианта- чистота

depth— относительная глубина, =2z/(x+y)

table— относительный размер вершины брилианта

price— цена

x— длина

y— ширина

z— высота

4. Типы признаков

X— порядковый

carat— количественные (дискретный признак, мода встречается 2604 раза)

cut— качественные

color— качественные

clarity— качественные

depth— количественные (дискретный признак, мода встречается 2239 раза)

table— количественные (дискретный признак, мода встречается 9881 раза)

price— количественные (ближе к непрерывному признаку, мода встречается 132 раза)

x— количественные (дискретный признак, мода встречается 448 раза)

y— количественные (дискретный признак, мода встречается 437 раза)

z— количественные (дискретный признак, мода встречается 767 раза)

summarize(df, across(carat:z, function(x) max(table(x))))
##   carat cut color clarity depth table price  x  y  z
## 1    55 394   187     247    46   187     5 12 12 18

5. Порядковые признаки

Здесь необходимо проверить соответствие текстовых меток порядкового признака (если такие использованы) их естественному порядку.
В данном датасете таких признаков нет.

6. Matrix plot, outliers, etc.

ggpairs(df, title="correlogram", columns=c(2,6:11), upper = list(continuous = "points"), diag = list(continuous = "barDiag"))

Наблюдаются outliers почти на всех скатерплотах. Удалим некторые слишком выделяющиеся значения, которые могут являться ошибками в данных.

dfo <- df
dfo[rownames(dfo)[dfo$x == 0 | dfo$z == 0 | dfo$table > 90 | dfo$z > 30 | dfo$y > 30 ],] <- NA
ggpairs(dfo, title="correlogram", columns=c(2,6:11), diag = list(continuous = "barDiag"))

7. Симметричность распределений.

Из matrixplot видно, что распределение carat, price, x, y, z - сильно несимметричное с хвостом вправо, поэтому прологарифмируем его и построим заново матрикс плот.

dfol <- transform(dfo, price=log(price), carat=log(carat), x=log(x), y=log(y), z=log(z))
names(dfol)[names(dfol) == 'price'] <- 'log_price'
names(dfol)[names(dfol) == 'carat'] <- 'log_carat'
names(dfol)[names(dfol) == 'x'] <- 'log_x'
names(dfol)[names(dfol) == 'y'] <- 'log_y'
names(dfol)[names(dfol) == 'z'] <- 'log_z'
ggpairs(dfol, title="correlogram", columns=c(2,6:11), diag = list(continuous = "barDiag"))

Распределения стали симметричнее, а зависимости более линейными.

8. Аутлаеры

Наибольшая линейная зависимость наблюдается между carat, x; carat, y; carat, z; x, y; x, z; y, z. Наименьшая- price, depth; x, depth; y, depth; z, depth. Удалим значения, которые могут быть особыми индивидами.

dfolo <- na.omit(dfol)
dfolo[rownames(dfolo)[dfolo$table > 65 | dfolo$table < 50 | dfolo$depth > 65 | dfolo$depth < 60 | dfolo$log_z < 0.5],] <- NA
ggpairs(dfolo, title="correlogram", columns=c(2,6:11), diag = list(continuous = "barDiag"))

10. Раскраска по категориальным признакам

Cut

ggpairs(dfolo, columns=c(2,6:8), ggplot2::aes(colour=cut), diag = list(continuous = "barDiag"))

Color

ggpairs(dfolo, columns=c(2,6:8), ggplot2::aes(colour=color), diag = list(continuous = "barDiag"))

Clarity

ggpairs(dfolo, columns=c(2,6:8), ggplot2::aes(colour=clarity), diag = list(continuous = "barDiag"))

Не было выявлено неоднородностей.

11. Descriptive statistics

summary(dfo[-1])
##      carat               cut          color        clarity        depth      
##  Min.   :0.2300   Fair     : 35   E      :187   SI1    :247   Min.   :56.30  
##  1st Qu.:0.3900   Good     :105   G      :185   SI2    :195   1st Qu.:61.10  
##  Median :0.7100   Ideal    :394   H      :170   VS2    :195   Median :61.90  
##  Mean   :0.8014   Premium  :253   F      :169   VS1    :150   Mean   :61.84  
##  3rd Qu.:1.0500   Very Good:212   D      :130   VVS2   :104   3rd Qu.:62.50  
##  Max.   :3.0100   NA's     :  1   (Other):158   (Other):108   Max.   :71.60  
##  NA's   :1                        NA's   :  1   NA's   :  1   NA's   :1      
##      table           price               x               y        
##  Min.   :51.00   Min.   :  390.0   Min.   :3.920   Min.   :3.950  
##  1st Qu.:56.00   1st Qu.:  923.5   1st Qu.:4.690   1st Qu.:4.690  
##  Median :57.00   Median : 2655.0   Median :5.720   Median :5.740  
##  Mean   :57.53   Mean   : 3862.6   Mean   :5.733   Mean   :5.735  
##  3rd Qu.:59.00   3rd Qu.: 5223.5   3rd Qu.:6.555   3rd Qu.:6.570  
##  Max.   :70.00   Max.   :18717.0   Max.   :8.990   Max.   :8.930  
##  NA's   :1       NA's   :1         NA's   :1       NA's   :1      
##        z        
##  Min.   :1.070  
##  1st Qu.:2.900  
##  Median :3.550  
##  Mean   :3.545  
##  3rd Qu.:4.050  
##  Max.   :5.860  
##  NA's   :1
summary(na.omit(dfolo[-1]))
##    log_carat               cut      color      clarity        depth      
##  Min.   :-1.46968   Fair     : 10   D:118   SI1    :223   Min.   :60.00  
##  1st Qu.:-0.96758   Good     : 91   E:167   VS2    :176   1st Qu.:61.30  
##  Median :-0.35667   Ideal    :386   F:145   SI2    :166   Median :62.00  
##  Mean   :-0.41478   Premium  :210   G:168   VS1    :135   Mean   :61.98  
##  3rd Qu.: 0.04401   Very Good:194   H:153   VVS2   : 98   3rd Qu.:62.50  
##  Max.   : 0.94391                   I: 95   VVS1   : 54   Max.   :65.00  
##                                     J: 45   (Other): 39                  
##      table         log_price         log_x           log_y      
##  Min.   :53.00   Min.   :5.969   Min.   :1.366   Min.   :1.374  
##  1st Qu.:56.00   1st Qu.:6.810   1st Qu.:1.538   1st Qu.:1.539  
##  Median :57.00   Median :7.822   Median :1.739   Median :1.740  
##  Mean   :57.32   Mean   :7.754   Mean   :1.719   Mean   :1.720  
##  3rd Qu.:59.00   3rd Qu.:8.543   3rd Qu.:1.876   3rd Qu.:1.876  
##  Max.   :65.00   Max.   :9.837   Max.   :2.177   Max.   :2.167  
##                                                                 
##      log_z       
##  Min.   :0.8671  
##  1st Qu.:1.0578  
##  Median :1.2613  
##  Mean   :1.2409  
##  3rd Qu.:1.3962  
##  Max.   :1.6827  
## 
summarize(na.omit(dfolo[-1]), across(c(log_carat, depth:log_z), list(kurtosis = kurtosis, skewness = skewness)))
##   log_carat_kurtosis log_carat_skewness depth_kurtosis depth_skewness
## 1           1.834324         0.05581454       3.033115      0.3349012
##   table_kurtosis table_skewness log_price_kurtosis log_price_skewness
## 1       3.595072      0.6033307           1.868285          0.1150333
##   log_x_kurtosis log_x_skewness log_y_kurtosis log_y_skewness log_z_kurtosis
## 1       1.864098     0.06898368       1.850999     0.06982123       1.837582
##   log_z_skewness
## 1      0.0561814

kurtosis и skewness не равен 0 ни у одного признака, можно предположить, что выборки не из нормального распределения, позже это будет проверено при помощи теста Шапиро-Уилка.

У прологорифмированных признаков медиана и математическое ожидание стали ближе. Квартили, а также min, max графически изображены при помощи Boxplot и рассмотрены в пункте 2.2.

2.1 Выбор категоризующей переменной

В качестве категоризующего признака возьмем cut. Сравнивать будем Ideal и Premium, как имеющие наибольшее количество наблюдений.

library(ppcor) #Библиотека, позволяющая проверять гипотезы о значимости коэффициента частной корреляции.
library(Hmisc) #Библиотека, позволяющая строить корреляционные матрицы вместе с матрицами соответствующих p-value теста о значимости коэффициента корреляции.
library(ggpubr)

2.2 Boxplot

dfcomp <- dfolo %>% filter(cut == "Ideal" | cut == "Premium")
dfcomp$cut <- droplevels(dfcomp$cut)
bwplot(log_carat ~ cut, data = dfcomp, col = c("forestgreen", "gold"), main = "log_carat", xlab = "cut")

bwplot(depth ~ cut, data = dfcomp, col = c("forestgreen", "gold"), main = "depth", xlab = "cut")

bwplot(table ~ cut, data = dfcomp, col = c("forestgreen", "gold"), main = "table", xlab = "cut")

На boxplot признака log_carat: Примерно равный разброс, но у Premium медиана больше.

На boxplot признака depth: Примерно равные медианы, но у Premium больше разброс, наблюдается несколько outliers.

На boxplot признака table: Примерно равный разброс, но у Premium медиана больше, наблюдается несколько outliers.

2.3 Нормальность признаков

\(\left\{ \left(x_{i},cdf_{0}^{-1}\left(\widehat{cdf}_{n}(x_{i})+\frac{1}{2n}\right)\right)\right\} _{i=1}^{n}.\)

Частный случай Q-Q plot для \(cdf_{0}^{-1}=cdf_{N(0,1)}^{-1}\) называется normal probability plot.

Если \(\hat{P}_{n}->P_{\xi}\), то оба графика будут стремиться к \(y=x\). Референсной прямой normal probability plot будет \(y=\sqrt{\widehat{D\xi}}\cdot x+\widehat{E\xi}\).

Тест Шапиро-Уилка проверяет нулевую гипотезу о том, что выборка из нормального распределения

ggqqplot(subset(dfcomp, cut == "Ideal")$log_carat, ylab = "log_carat")

ggqqplot(subset(dfcomp, cut == "Ideal")$depth, ylab = "depth")

ggqqplot(subset(dfcomp, cut == "Ideal")$table, ylab = "table")

shapiro.test(subset(dfcomp, cut == "Ideal")$log_carat)
## 
##  Shapiro-Wilk normality test
## 
## data:  subset(dfcomp, cut == "Ideal")$log_carat
## W = 0.93548, p-value = 6.958e-12
shapiro.test(subset(dfcomp, cut == "Ideal")$depth)
## 
##  Shapiro-Wilk normality test
## 
## data:  subset(dfcomp, cut == "Ideal")$depth
## W = 0.97893, p-value = 2.068e-05
shapiro.test(subset(dfcomp, cut == "Ideal")$table)
## 
##  Shapiro-Wilk normality test
## 
## data:  subset(dfcomp, cut == "Ideal")$table
## W = 0.93437, p-value = 5.25e-12
#Далее, аналогично для другого качества огранки

ggqqplot(subset(dfcomp, cut == "Premium")$log_carat, ylab = "log_carat")

ggqqplot(subset(dfcomp, cut == "Premium")$depth, ylab = "depth")

ggqqplot(subset(dfcomp, cut == "Premium")$table, ylab = "table")

shapiro.test(subset(dfcomp, cut == "Premium")$log_carat)
## 
##  Shapiro-Wilk normality test
## 
## data:  subset(dfcomp, cut == "Premium")$log_carat
## W = 0.93599, p-value = 5.702e-08
shapiro.test(subset(dfcomp, cut == "Premium")$depth)
## 
##  Shapiro-Wilk normality test
## 
## data:  subset(dfcomp, cut == "Premium")$depth
## W = 0.95925, p-value = 1.016e-05
shapiro.test(subset(dfcomp, cut == "Premium")$table)
## 
##  Shapiro-Wilk normality test
## 
## data:  subset(dfcomp, cut == "Premium")$table
## W = 0.92838, p-value = 1.336e-08

Гипотеза о нормальности распределения отвергается как при оценке p-value, так и при рассмотрение Q-Q plot. Так же на Q-Q plot можно заметить дискретность распределений, особенно у table.

2.4 t-test, критерий Манна-Уитни

Двухвыборочный \(t\)-критерий:

\(H_{0}:E\xi_{1}=E\xi_{2}\).

\(t=\frac{\bar{x}-\bar{y}}{\sqrt{D(\bar{x}-\bar{y})}}\xrightarrow{\sim}N(0,1).\)

t-test более мощный против гипотезы о разных математических ожиданиях, а тест wilcox имеет большую устойчивость, так как является ранговым тестом и он проигнорирует выбросы.

t.test(log_carat ~ cut, data = dfcomp)
## 
##  Welch Two Sample t-test
## 
## data:  log_carat by cut
## t = -6.4454, df = 396.34, p-value = 3.353e-10
## alternative hypothesis: true difference in means between group Ideal and group Premium is not equal to 0
## 95 percent confidence interval:
##  -0.4220803 -0.2247763
## sample estimates:
##   mean in group Ideal mean in group Premium 
##            -0.5803133            -0.2568850
t.test(depth ~ cut, data = dfcomp)
## 
##  Welch Two Sample t-test
## 
## data:  depth by cut
## t = 1.2824, df = 345.41, p-value = 0.2006
## alternative hypothesis: true difference in means between group Ideal and group Premium is not equal to 0
## 95 percent confidence interval:
##  -0.04338794  0.20596381
## sample estimates:
##   mean in group Ideal mean in group Premium 
##              61.75415              61.67286
t.test(table ~ cut, data = dfcomp)
## 
##  Welch Two Sample t-test
## 
## data:  table by cut
## t = -23.527, df = 384.31, p-value < 2.2e-16
## alternative hypothesis: true difference in means between group Ideal and group Premium is not equal to 0
## 95 percent confidence interval:
##  -2.968272 -2.510410
## sample estimates:
##   mean in group Ideal mean in group Premium 
##              55.94637              58.68571
wilcox.test(log_carat ~ cut, data = dfcomp)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  log_carat by cut
## W = 28214, p-value = 8.531e-10
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(depth ~ cut, data = dfcomp)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  depth by cut
## W = 42160, p-value = 0.4167
## alternative hypothesis: true location shift is not equal to 0
wilcox.test(table ~ cut, data = dfcomp)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  table by cut
## W = 5975.5, p-value < 2.2e-16
## alternative hypothesis: true location shift is not equal to 0

Гипотеза о равенстве средних отвергается для всех подвыборок по cut.

2.5 Критерий Колмогорова-Смирнова

Рассматривается \(H_{0}:P_{\xi_{1}}=P_{\xi_{2}}\) против \(H_{1}:P_{\xi_{1}}\neq P_{\xi_{2}}\) и оба распределения абсолютно непрерывны. В качестве статистики используется \(D=\sup_{x}\left|\widehat{cdf}_{\xi_{1}}(x)-\widehat{cdf}_{\xi_{2}}(x)\right|.\)

Критерий Колмогорова-Смирнова применим для непрерывных признаков, поэтому исключим из рассмотрения table

ks.test(dfcomp[dfcomp$cut == "Ideal", 2], dfcomp[dfcomp$cut == "Premium", 2]) 
## 
##  Two-sample Kolmogorov-Smirnov test
## 
## data:  dfcomp[dfcomp$cut == "Ideal", 2] and dfcomp[dfcomp$cut == "Premium", 2]
## D = 0.29507, p-value = 1.037e-10
## alternative hypothesis: two-sided
ks.test(dfcomp[dfcomp$cut == "Ideal", 6], dfcomp[dfcomp$cut == "Premium", 6]) 
## 
##  Two-sample Kolmogorov-Smirnov test
## 
## data:  dfcomp[dfcomp$cut == "Ideal", 6] and dfcomp[dfcomp$cut == "Premium", 6]
## D = 0.13161, p-value = 0.01799
## alternative hypothesis: two-sided

Гипотезы о равенстве распределений отвергаются.

3.1 Анализ зависимостей по группе.

Посмотрим на матрикс плот данных.

Cut

ggpairs(dfolo, columns=c(2,6:8), ggplot2::aes(colour=cut), diag = list(continuous = "barDiag"))

Color

ggpairs(dfolo, columns=c(2,6:8), ggplot2::aes(colour=color), diag = list(continuous = "barDiag"))

Clarity

ggpairs(dfolo, columns=c(2,6:8), ggplot2::aes(colour=clarity), diag = list(continuous = "barDiag"))

Корреляция Пирсона измеряет линейную зависимость между двумя переменными (x и y).

Мера линейной зависимости между случайным величинами \(\xi\) и \(\eta\) есть коэффициент корреляции Пирсона \(\rho=\frac{cov(\xi,\eta)}{\sqrt{D\xi}\sqrt{D\eta}}.\)

3.2 Коэффициент корреляции Пирсона.

cor.test( ~ log_carat + depth, data = dfolo, subset = cut == "Ideal", method = "pearson") #Считает коэффициент корреляции Спирмена между двумя столбцами матрицы.
## 
##  Pearson's product-moment correlation
## 
## data:  log_carat and depth
## t = 1.2066, df = 384, p-value = 0.2283
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.03859601  0.16028946
## sample estimates:
##        cor 
## 0.06145674
cor.test( ~ depth + table, data = dfolo, subset = cut == "Ideal", method = "pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  depth and table
## t = -5.2785, df = 384, p-value = 2.183e-07
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.3508070 -0.1645548
## sample estimates:
##        cor 
## -0.2600986
cor.test( ~ table + log_carat, data = dfolo, subset = cut == "Ideal", method = "pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  table and log_carat
## t = 2.5335, df = 384, p-value = 0.01169
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.02877438 0.22515639
## sample estimates:
##       cor 
## 0.1282221
rcorr(as.matrix(subset(dfolo, cut == 'Ideal', select=c(log_carat, depth, table))), type = "pearson")
##           log_carat depth table
## log_carat      1.00  0.06  0.13
## depth          0.06  1.00 -0.26
## table          0.13 -0.26  1.00
## 
## n= 386 
## 
## 
## P
##           log_carat depth  table 
## log_carat           0.2283 0.0117
## depth     0.2283           0.0000
## table     0.0117    0.0000
cor.test( ~ log_carat + depth, data = dfolo, subset = cut == "Premium", method = "pearson") #Считает коэффициент корреляции Пирсона между двумя столбцами матрицы.
## 
##  Pearson's product-moment correlation
## 
## data:  log_carat and depth
## t = 0.094237, df = 208, p-value = 0.925
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.1289705  0.1417990
## sample estimates:
##         cor 
## 0.006534022
cor.test( ~ depth + table, data = dfolo, subset = cut == "Premium", method = "pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  depth and table
## t = -0.26166, df = 208, p-value = 0.7938
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.1531539  0.1175395
## sample estimates:
##         cor 
## -0.01813964
cor.test( ~ table + log_carat, data = dfolo, subset = cut == "Premium", method = "pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  table and log_carat
## t = 0.6197, df = 208, p-value = 0.5361
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.09300244  0.17728855
## sample estimates:
##        cor 
## 0.04292854
rcorr(as.matrix(subset(dfolo, cut == 'Premium', select=c(log_carat, depth, table))), type = "pearson")
##           log_carat depth table
## log_carat      1.00  0.01  0.04
## depth          0.01  1.00 -0.02
## table          0.04 -0.02  1.00
## 
## n= 210 
## 
## 
## P
##           log_carat depth  table 
## log_carat           0.9250 0.5361
## depth     0.9250           0.7938
## table     0.5361    0.7938

Не отвергается гипотеза о том, что корреляция = 0 между log_carat and depth для подвыборки cut == “Ideal”, log_carat and depth для подвыборки cut == “Premium” и table and log_carat для подвыборки cut == “Premium”. Значения коэффициента корреляции приведены в таблицах.

3.3 Коэффициент корреляции Спирмена

Выборочный коэффициент Спирмана: \(\hat{\rho}_{S}=\frac{1/n\cdot\sum_{i=1}^{n}R_{i}T_{i}-\bar{R}\bar{T}}{\sqrt{1/n\cdot\sum_{i=1}^{n}\left(R_{i}-\bar{R}\right)^{2}}\sqrt{1/n\cdot\sum_{i=1}^{n}\left(T_{i}-\bar{T}\right)^{2}}}.\)

Если нет повторяющихся наблюдений, то знаменатель будет одним и тем же у всех выборок объема \(n\), значит его можно посчитать заранее. В этом (и только этом) случае, справедлива более простая формула: \(\hat{\rho}_{S}=1-\frac{6\sum_{i=1}^{n}(R_{i}-T_{i})^{2}}{n^{3}-n}.\)

Коэффициент корреляции Спирмена являеется ранговым и соответсвенно устойчивым к выбросам. В случае если распределение нормальное коэффициент Спирмена и коэффициент Пирсона измеряют одно и то же.

cor.test( ~ log_carat + depth, data = dfolo, subset = cut == "Ideal", method = "spearman") #Считает коэффициент корреляции Спирмена между двумя столбцами матрицы.
## 
##  Spearman's rank correlation rho
## 
## data:  log_carat and depth
## S = 8831945, p-value = 0.1232
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##        rho 
## 0.07859913
cor.test( ~ depth + table, data = dfolo, subset = cut == "Ideal", method = "spearman")
## 
##  Spearman's rank correlation rho
## 
## data:  depth and table
## S = 11718570, p-value = 1.016e-05
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##        rho 
## -0.2225507
cor.test( ~ table + log_carat, data = dfolo, subset = cut == "Ideal", method = "spearman")
## 
##  Spearman's rank correlation rho
## 
## data:  table and log_carat
## S = 8764740, p-value = 0.09303
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##        rho 
## 0.08561034
rcorr(as.matrix(subset(dfolo, cut == 'Ideal', select=c(log_carat, depth, table))), type = "spearman")
##           log_carat depth table
## log_carat      1.00  0.08  0.09
## depth          0.08  1.00 -0.22
## table          0.09 -0.22  1.00
## 
## n= 386 
## 
## 
## P
##           log_carat depth  table 
## log_carat           0.1232 0.0930
## depth     0.1232           0.0000
## table     0.0930    0.0000
cor.test( ~ log_carat + depth, data = dfolo, subset = cut == "Premium", method = "spearman") #Считает коэффициент корреляции Спирмена между двумя столбцами матрицы.
## 
##  Spearman's rank correlation rho
## 
## data:  log_carat and depth
## S = 1526959, p-value = 0.8776
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##        rho 
## 0.01069433
cor.test( ~ depth + table, data = dfolo, subset = cut == "Premium", method = "spearman")
## 
##  Spearman's rank correlation rho
## 
## data:  depth and table
## S = 1581454, p-value = 0.7229
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##         rho 
## -0.02461259
cor.test( ~ table + log_carat, data = dfolo, subset = cut == "Premium", method = "spearman")
## 
##  Spearman's rank correlation rho
## 
## data:  table and log_carat
## S = 1422979, p-value = 0.2601
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##       rho 
## 0.0780617
rcorr(as.matrix(subset(dfolo, cut == 'Premium', select=c(log_carat, depth, table))), type = "spearman")
##           log_carat depth table
## log_carat      1.00  0.01  0.08
## depth          0.01  1.00 -0.02
## table          0.08 -0.02  1.00
## 
## n= 210 
## 
## 
## P
##           log_carat depth  table 
## log_carat           0.8776 0.2601
## depth     0.8776           0.7229
## table     0.2601    0.7229

Не отвергается гипотеза о том, что корреляция = 0 между log_carat and depth для подвыборки cut == “Ideal”, log_carat and depth для подвыборки cut == “Premium” и table and log_carat для подвыборки cut == “Premium”. Значения коэффициента корреляции приведены в таблицах.

Предложим, что вес в каратах является внешним фактором, влияющим на корреляцию между table и depth. Для проверки этого предположения посмотрим на коэффициент частной корреляции.

3.4 Причины, следствия и частные корреляции.

cor.test( ~ table + depth, data = dfolo, subset = cut == "Ideal", method = "pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  table and depth
## t = -5.2785, df = 384, p-value = 2.183e-07
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.3508070 -0.1645548
## sample estimates:
##        cor 
## -0.2600986
cor.test( ~ table + depth, data = dfolo, subset = cut == "Premium", method = "pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  table and depth
## t = -0.26166, df = 208, p-value = 0.7938
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.1531539  0.1175395
## sample estimates:
##         cor 
## -0.01813964
dft <- subset(dfolo, cut == "Ideal")
pcor.test(dft$depth, dft$table, dft$log_carat, method = "pearson") #Смотрим на проверку гипотезы о значимости коэффициента частной корреляции
##     estimate      p.value statistic   n gp  Method
## 1 -0.2707209 6.824109e-08  -5.50363 386  1 pearson
dfpart <- na.omit(subset(dfolo, cut == "Premium"))
pcor.test(dfpart$depth, dfpart$table, dfpart$log_carat, method = "pearson")
##      estimate   p.value  statistic   n gp  Method
## 1 -0.01843753 0.7910309 -0.2653149 210  1 pearson

Значение коэффициента частной корреляции слабо отличается от исходного коэффициента, можно сделать вывод, что корреляция между table и depth обеспечена внутренними связями между признаками, это поддтверждается и значениями p-value.